Stochastic Calculus for Assets with Non-Gaussian Price Fluctuations * 
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From the path integral formalism for price fluctuations with non-Gaussian distribu- 
tions I derive the appropriate stochastic calculus replacing Ito's calculus for stochastic 
fluctuations. 
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The logarithms of assets prices in financial markets do not fluctuate with Gaussian distributions. This 
was first noted by Pareto in the 19th century |l| and emphasized by Mandelbrot in the 1960s 0]- The 
true distributions possess large tails which may be approximated by various other distributions, most 
prominently the truncated Levy distributions [p|- p|, the Meixner distributions H0, or the generalized 
hyperbolic distributions and their descendents |8H32| . The stochastic differential equations associated 
with such distributions cannot be treated with the popular Ito calculus. The purpose of this paper 
is to develop an appropriate calculus to replace it. In Section Owe briefly recapitulate the Gaussian 
approximation and set the stage for the generalization in Section [III]. 
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II. GAUSSIAN APPROXIMATION TO FLUCTUATION PROPERTIES OF STOCK PRICES 

Let S(t) denote the price of some stock. Over long time spans, the average over many stock prices 
has a time behavior which can be approximated by pieces of exponentials. This is why they are usually 
plotted on a logarithmic scale. For an illustration see the Dow- Jones industrial index over 60 years in 
Fig. | 
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FIG. 1. Periods of exponential growth of price index averaged over major industrial stocks in the Unites States 
over 60 years. 



For a liquid market with many participants, the price fluctuations seem to be driven by a stochastic noise 
with a white spectrum, as illustrated in Fig. @. 
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FIG. 2. Fluctuation spectrum of exchange rate DM/US$ as a function of the frequency, indicating a white 
noise in the stochastic differential equation (Q) (from the textbook [^]). 

The fluctuations of the index have a certain average width called the volatility of the market. Over 
longer time spans, the volatility changes stochastically, as illustrated by the data of the S&P500 index 
over the years 1984-1997 shown in Fig. ||. In particular, there are strong increases short before market 
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FIG. 3. (a) The S&P 500 index for the 13-year period 1 Jan 1984 - 31 Dec 1996 at interval of 1 min. Large 
fluctuations appeared on 19 Oct 1987 (black Monday), (b) Volatility over times T=lmon (8190 min) for courses 
taken every 30min. The precursors of the '87 crash are indicated by arrows (from Ref. J3^|). 

crashes. The theory to be developed here will ignore these fluctuations and assume a constant volatility. 
For recent work taking them into account see [M. 

The distribution of the logarithms of the volatilities is normal as shown in Fig. || 

Probability xl0~ 3 

1.0 
0.8 
0.6 
0.4 
0.2 



0.0 
0.000 



r \ 


\ + -300-min-data 


ji \ 


\ — dog-normal-fit 


f 


\ \ — -Gaussian-fit 

-K \ 
\\ 
4v\ 




V 

\ 
A 




\ \ 
\ A 


J 


\ A^ 

,"^!^A^ 







0.003 



0.004 



0.001 0.002 
Votatility a 

FIG. 4. Comparison of the best log-normal and Gaussian fits to volatilities over 300-min (from Ref. ]33|). 
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An individual stock will in general have larger volatility than an average market index, especially when 
the company is small and only few stocks are traded per day. 

To lowest approximation, the stock price S(t) satisfies the simplest stochastic differential equation for 
exponential growth 



S(t) 



where 77(f) is a white noise with the correlation functions 

( V (t)) =0, ( V (t)r ] (t')) = a 2 S(t-t'). 



(1) 



(2) 



The standard deviation a parametrizes the volatility of the stock price, which is measured by the expec- 
tation value 





~S(t)~ 


( 


S(t) 



-I 2 



dt = a 2 . 



The logarithm of the stock price 



x(t) = log S(t) 



(3) 



(4) 



does not simply satisfy the stochastic linear-growth differential equation x(t) = S(t)/S(t) = r$ + Tf{t), 
but possesses another growth rate r x : 



x(t) = r x + r)(t). 



(5) 



The change of the growth rate is due the stochastic nature of x(t) and S(t). According to Ito's rule, a 
function of a stochastic variable satisfies the relation 



/(*(*)) = d x f(x(t))x(t) + — f" (x(t)). 



(6) 



This is derived by expanding f(x(t + e)) in a Taylor series 



f(x(t + e)) = f(x(t)) + f'(x(t)) / dt'x(t') 



l+c 



+ ^/ (3) (^W) 



l+e- 



dti 

(+< 



t+e 



dt 2 x(ti)x(t 2 ) 



(7) 



dh 



L + c 



dh 



t+c- 



tzx{ti)x(t 2 )x(ts,) + 



replacing the higher correlation functions on the right-hand side by their expectation values, and observing 
that only the quadratic noise term in x(t) = r x + rj(t) contributes to linear order in e. For the logarithmic 
function (Q), the expansion yields 



S(t) 2 



S(t) 



dt + 



(8) 



which becomes, after the replacement of the quadratic term by its expectation value (||), 

./ \ S(t) 1 2 

Inserting here Eqs. (0), and (||), we find that the linear growth rate of x(t) is related to the exponential 
growth rate of S{t) by 



(9) 
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r x =r s -^a 2 . (10) 

In praxis, this relation implies that if we fit a straight line through a plot of the logarithms of stock prices, 
the forward extrapolation of the average stock price is given by 

(S{t)} = S(0) e rst = S{0) e (^+ CT2 /2)t. ( U ) 

A typical set of solutions of the stochastic differential equation (||) is shown in Fig. ||. 




FIG. 5. Behavior of logarithm of stock price following the stochastic differential equation (hi). 



III. LEVY DISTRIBUTIONS 



The description of the fluctuations of the logarithms of the stock prices around the linear trend by a 
Gaussian distribution is only a rough approximation to the real stock prices. As explained before, these 
have volatilities depending on time. More severely, if observed in smaller time intervals, for instance 
every day or hour, they have distributions in which rare events have a much higher probability than in 
Gaussian distributions, whose exponential tails are extremely small. Following Pareto JjJ, Mandelbrot 
emphasized in the 1960s that fluctuations of financial assets could be fitted much better with the help of 
Levy distributions M . These distributions are defined by the Fourier transform 



(12) 



with 



(13) 



The Gaussian distribution are recovered in the limit p, — > 2. 

For large x, the Levy distribution (|l2|) falls off with the characteristic power behavior 



(14) 



These power falloffs are referred to as Paretian tails of the distributions. The amplitude of the tails is 
found by approximating the integral (12) for large x, where only small momenta contribute, as follows 



dp 
27 



At, 



(15) 



with 



2p, 



dp' 



p'V" cos p' 



27T/1 



sin(7r/i/2)r(l + M)- 



(16) 



The stock market data are fitted best with p, between 1.2 and 1.5 ||, and we shall use p, = 3/2 most 
of the time for simplicity, where one has 



,3/2 _ 1 

^ "4V2V 



(17) 
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The full Taylor expansion of the Fourier transform ( J13| ) yields the asymptotic series 

n=0 17 u n=0 1 1 

This series is not useful for practical calculations since it fails to reproduce the exponential tails of the 
distribution, a typical shortcoming of large- a; expansions. In particular, it does not reduce to the Gaussian 
distribution in the limit ji — ► 2. 



A. Truncated Levy Distributions 

An undesirable property of the Levy distributions which is incompatible with financial data is that 
their fluctuation width diverges for /i < 2, since 



Jj XX 2 L» 2 { X )=-—L^{p) 



a 2 = (x 2 )= I J~~ 2 TM 



2 

(19) 

p=0 



is infinite. In contrast, real stock prices have a finite width. To account for both Paretian tails and finite 
width one may introduce the so-called truncated Levy distributions ||. They are defined by 

L<$ a \x)= £^&r\v)^\ (20) 

with a Fourier transform L 2 (p) which generalizes the function (|l^). It is conveniently written as an 
exponential of a "Hamiltonian function" H(p): 

L$r>(p) = e- H M, (21) 

with 

a 2 a 2 ~» 



H(p) = — — r [(a + ipY + (a~ ipY - 2a"} 

2 W 

2 (a 2 +p 2 ) >1 / 2 cos[n arctan(p/a)] — a" 



(22) 



The asymptotic behavior of the truncated Levy distributions differs from the power behavior of the 
Levy distribution in Eq. ( |l5| ) by an exponential factor e~ ax , which guarantees the finiteness of the width a 
and of all higher moments. An estimate of the leading term is again obtained from the Fourier transform 

L[T\x) « e 2Ca " [°° ^{l + CKa + ipr + ia-ipr^e^ 

J -co 27T 

^rd+M)^^, (23) 



x— >oo 



where 



Cs tmwT (24) 

The integral follows directly from the formulas [H 

,2^ V ^ r(- M )xi+^' 7_ 0O 2vr V w r(- M ) |a:|i+M' 1 j 

and the identity for Gamma functions l/T(—z) = —T(l + z) sm(ir z) / it . The full expansion is integrated 
with the help of the formula |3q] 
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/ $-{a + i P Y{a-ip) v e v 
J-oo 2tt 



-W / (^- Al )/2,(i+ A1 + 1 ,)/2(2aa;) x > 0, 



= i^r^ lx ' /2 V ^ (26) 

[ T(-v) ^-^72,(1+^+^/2 (2ag) x<0, 

where the Whittaker functions W / ( !y _ Al )/ 2 ,(i+M+i / )/2(2 a;c ) can be expressed in terms of Rummer's confluent 
hypergeometric function \Fi(a;b;z) as 

= r(yt-K-x) zK+1/2e ~ z/2lFl{1/2 + K ~ A; 2K + l]z) 

n2, '' ! -«-"+^V*/Vi(1/2-k-A;-2k+1;z). (27) 



r(i/2 + K-A)' 



For z/ = 0, only x > gives a nonzero integral (|2"6|), which reduces, with W_ AI / 2 ,i/2+^/2( z ) = z z ^ 2 : 
to the left equation in (|25|). Setting /x = v we find 



(28) 



Inserting 



Wo,i/2+ v («) = \—K 1/2+u {z/2), (29) 



7T 

we may write 

| (o . + ^ . (^)' /J+ "_l_ A - 1/2+ „ (tlW) ^ (30) 

For v = — 1 where i4T„ 1 / 2 (^) = ^i/2( z ) = y/'x/'2ze~ z , this reduces to 

^ l- e *« = -Le-«W. (31) 
27T or + #r 2a 

Summing up all terms in the expansion 

L<$*\x) « e 2C ^ |l + fj £g£ [(a + ip)" + (a - i P )*A e*» (32) 

yields the true asymptotic behavior which differs from the estimate ( p3| ) by a constant factor | f37f : 

- e ( 2 -^^r(i + M )^Ma^. (33) 

In contrast to Gaussian distributions which are characterized completely by the width a, the truncated 
Levy distributions contain three parameters er, fi, and a. Best fits to two types of fluctuating market 
prices are shown in Fig. ^, in which we plot the cumulative probabilities 

/Sx poo 
dx L%' a] (x) , P> (fa) = dx L [ £ a) (x) = 1 - P< (Sx) . (34) 
-co J Sx 

For negative price fluctuations 8x, the plot shows P < (5x), for positive price fluctuations P > (5x). By 
definition, P<(-oo) = 0, P<(0) = 1/2, P<(oo) = 1, and P>(-oo) = 1, P>(0) = 1/2, P>(oo) = 0. To fit 
the general shape, one chooses an appropriate parameter /i which turns out to be rather universal, close 
to /z = 3/2. The remaining two parameters fix all expansion coefficients of Hamiltonian (E3): 
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Hip) = \c2 p 2 - ^c 4 / + ^c 6 p 6 - ic 8 p s + . . . . (35) 
The numbers c n are referred to as the cumulants of the truncated Levy distribution. They are equal to 

C2 = c 2 , 

C4 = <t 2 (2 — //)(3 — [i)a~ 2 , 

c 6 = a 2 (2 - M )(3 - M )(4 - M )(5 - p)^ 4 , 

2 T(2n-/x) 2 _ 2 „ 



see extra fig6a.pdf see extra fig6b.pdf 



FIG. 6. Best fit of cumulative versions (|34|) of truncated Levy distribution, P< (Sx) for negative price fluctuations 
5x, and P>(Sx) for positive price fluctuations to the instantaneous fluctuations of the S&P 500 (l.h.s) and the 
DM/US$ exchange rate measured every fifteen minutes. The negative fluctuations lie on a slightly higher curve 
than the positive ones. The difference is often neglected. The parameters A and a are the size and truncation 
parameters of the distribution. The best value of /x is 3/2 (from |H|). 



The first determines the quadratic fluctuation width 

= c 2 =ct 2 , (37) 



p=0 



the second the fourth-order expectation 



= c 4 + 3c|, (38) 

p=0 



and so on: 

(x 6 ) = c 6 + 15c 4 c 2 + 154, (x s ) = c 8 + 28c 6 c 2 + 35c 2 , + 210c 4 c^ + 105(4 .... (39) 

In analyzing the data, one usually defines the so-called kurtosis, which is the normalized fourth-order 
cumulant 

K = c 4 = % = - 3. (40) 
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It depends on the parameters a, fi, a as follows 

_ (2-/i)(3-/i) 

K — 5 5 . 

<7 z or 

Given the volatility a and the kurtosis k, we extract the Levy parameter a from the equation 



(41) 



1/(2-^(3-,,) 
a V k 

In terms of k and a 2 , the expansion coefficients are 

_ _ 2 (5-^)(4- M ) _ _ 2 (7-//)(6-/x)(5-/j)(4-//) 
C4 - K ' C6 " K (3-,i)(2- M )' C8 " K (3-/*)>(2-/0» ' 



- = „n/2-l r (" ~ ^)/ r ( 4 - M) 

(3-M)«/2-2(2- M )»/2-2- 

For jit = 3/2, the second equation in ( (42| ) becomes simply 



(43) 



and the coefficients ( ^5[) 



^'^2 - c T 1 1 2 
C4 = K, C6 = — — K , Cg = 5 • 7 • 11 K , 

O 



= r(n-3/2)/r(5/2) /2 _! 

3«/2-2/2™-4 



(45) 



At zero kurtosis, the truncated Levy distribution reduces to a Gaussian distribution of width a. The 
change in shape for a fixed width and increasing kurtosis is shown in Fig. 




FIG. 7. Change in shape of truncated Levy distributions of width a = 1 with increasing kurtoses k = 
(Gaussian, solid curve), 1, 2,5, 10. 

From the S&P and DM/US$ data with time intervals At = 15min one extracts a 2 — 0.280 and 0.0163, 
and the kurtoses k = 12.7 and 20.5, respectively. This implies M w 3.57, a w 0.46 and M w 61.35, 
a w 1.50, respectively. 

The other normalized cumulants (eg, cs, . . .) are then all determined to be (1881.72, 788627.46, . . .) and 
(—4902.92, 3.3168 x 10 6 ,...), respectively. The cumulants increase rapidly showing that the expansion 
needs resummation. 

From the data, the other normalized cumulants are found by evaluating the ratios of expectation values 
- - ^ 1^ +™ 

6 ~(.T 2 )3 iJ (. x 2)2+ JU > 

(x 8 ) (x 6 ) (x 4 ) 2 (x 4 ) , , 

c " 8 = W ~ 28 FF ~ 35 W + 420 W " 630 ' ■ ■ ' ■ (46) 
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B. Asymmetric Truncated Levy Distributions 



We have seen in the data of Fig. || that the price fluctuations have a slight asymmetry: Price drops 
are slightly larger than rises. This is accounted for by an asymmetric truncated Levy distribution. It has 
the general form M 



LS B| %) = e- fl W 1 



with a Hamiltonian function 



H{p) 



[(a + ip) x (l + + (a- ip) x (l -13)- 2a x 



a 2 a 2 ~ x 
2 A(l-A) 

2 (a 2 + p 2 ) x / 2 {cos[A arctan(p/a)] + i/3sin[A arctan(p/a)]} — a x 



a x - 2 \{l-X) 



This has a power series expansion 



H(p) = icip+ ^c 2 p 2 - ^ C 3P 3 - tt c 4P 4 + «7TC 5 p 5 



1 

if 



1 

5? 



which differs from (|35|) by the extra odd coefficients 

9 a 



Ci = a 



a, 



(l-A)' 
c 3 = a 2 (2-\)a~ 1 p, 
eg = a 2 {2 - A) (3 - A)(4 - A)a~ 3 /?, 



(47) 



(48) 



(49) 



(50) 



2 r(2n+l_A) , 2 
c 2 n +1 = o r(2 _ A) a ft 

The general formula valid for even and odd n is 

o r ( n - A ) 2-« / 1 , " = even 



r(2 - A) 

In addition to the even expectation values (| 



j3 n = odd. 



(51) 



(52) 



there are now also odd expectation values: 



dxxL y T' ,J >(x) = i^-e- H ^ 



d 

dp 



ci, 



p=0 



(x 2 ) 



dxx 2 L { \- a ^ 



„3\ — 



dxx 3 L (X 2 a ' \x)= -i-^-e 
a d A p 



-H(p) 



c 2 + cf, 



p=0 



-H(p) 



C.3 + 3C 2 C! + c\, 



p=0 



dxx 4 L (X ' a ^(x) 



d 4 p 



-H(p) 



— c 4 + 4c 3 ci + ic\ + 6c 2 c 2 + c l5 



p=0 



The distribution is now centered around a nonzero average value: 

fi={x)= c x . 

The fluctuation width is given by 

a 2 ^{x 2 )-(x) 2 = ((x-{x)) 2 )=c 2 . 



(53) 
(54) 

(55) 
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The asymptotic behavior of the asymmetric truncated Levy distributions is given by a straightforward 
modification of (25): 



(a) 



2tt S 



a 2 a 2 ~ A 



2 A(l-A) 



[(a + ^p) A (l+/3) + (a - zp) A (l-/3) - 2a A 



rr^oo 7T \X\ + 

In analyzing the data, one now uses the skewness 

-gj = ^ ((x 3 ) - 3c 2Cl - c?) = ^ ((x 3 ) - 3( (z - (x» 2 )<x) - (x} 3 ) 



(56) 



C3 



(57) 



The kurtosis reads now 

C4 



At = C4 



1 



= — ((x 4 ) - 4c 3 ci - 3c 2 . - 6c 2 c 2 - c\) 



1 



(x 4 ) - 4((x - <x)) 3 )(x) - 3 ((x 2 ) - (xfY 6 «x 2 ) - <x) 2 ) <x) 2 - <x>^ 
(x-(x)) 4 )-3( 3 



x 2 - (x) 2 



It depends on the parameters a, A, /3, and a or At as follows 

(2 - \)0 



era 



The kurtosis has still the same dependence (ttlj) on a, A, a as before, such that we may also write 



'(2- X)k 



(3. 



V (3 -A) 

The average position = (x) is proportional to the ratio of skewness and kurtosis: 

3- A s 



(x) = C\ = CT Z 



(1-A) 



(3 = a 



1-A/c 



(58) 



(59) 



(60) 



(61) 



From the data one extracts volatility a, kurtosis At, and skewness s which determine completely the 
asymmetric truncated Levy distribution. The data are then plotted against x — (x) = x — p,, so that they 
are centered at the average position. This centered distribution will be denoted by L7 a '^(i), i.e. 



L^\x)=L^\ X -p) 
The Hamiltonian associated with this zero-average distribution is 

H(p) = H[p) - H'(0)p, 



(62) 



(63) 



and its expansion in power of the momenta starts out with p 2 , i.e., the first term in ( [i"9"| ) is subtracted. 
In terms of At, tr, and s, the normalized expansion coefficients are 



l/a _ x r(n-A)/r(4-A) 
(3 - A) n / 2_2 (2 - A)™/ 2 



1 



for 



V(3-A)/(2-A)«; 



n = even, 



n = odd. 



(64) 



The change in shape of the distributions of a fixed width and kurtosis with increasing skewness is shown 
in Fig. @. We have plotted the distributions centered around the average position x = (ci) which means 
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that we have removed from H(p) in (47), (|48|), and ( f49| ) the linear term icip. This subtracted Hamiltonian 
whose power series expansion begins with the term c 2 p 2 /2 will be denoted by 



H(p) = H{p) - H'(0) = \c 2 p 2 - t^c 3P 3 



1 

4! 



CAP 



(65) 




FIG. 8. Change in shape of truncated Levy distributions of width a = 1 and kurtosis k = 1 with increasing 
skewness s — (solid curve), 0.4, 0.8 . The curves are centered around (x). 



C. Meixner Distributions 

Quite reasonable fits to financial data are provided by the Meixner distributions |^J^] . They have the 
advantage that they can be given explicitly in configuration and momentum space 

^ {x) = [ \^r(2i) 2d |r {d + lx,a) |2 exp [bx/a] ' (66) 

, N f cos(6/2) ) 2d 

They have the same tail behavior as the truncated Levy distributions 

Mix) -> C±\x\ p e- a±lxl for x->±oo, (68) 

with 



[2cOs(b/2)] 2d 2lT 2 ;T,/iani/, 2, « . , _ ,,,/ , (;(n 

C±= 2a.r(2d) a^-i 6 1 ^ J , P - 2d - 1, ^ = (tt ± 6)/a. (69) 



The moments are 



/i=odtan(6/2), cr 2 = a 2 d/2 cos 2 (6/2), 

s=-v / 2sin(6/2)/v / d, k = [2 - cosb]/d, (70) 

such that we can calculate the parameters from the moments as follows: 

a 2 (2k -3s 2 ), d= 1 2 , 6 = -2arcsin^sVW2) ■ (71) 



a 2 



The Meixner distributions can be fitted quite well to the truncated Levy distribution in the regime 
of large probability. In doing so we observe that the variance a 2 and the kurtosis k are not the best 
parameters to match the two distributions. The large-probability regime of the distributions can be 
matched perfectly by choosing, in the symmetric case, the value and the curvature at the origin to be 
the same in both curves. This is seen in Fig. |9|. 

In the asymmetric case we have to match also the first and third derivatives. The derivatives of the 
Meixner distribution are: 
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M(0) 



M'(0) = 6 
M"(0) = - 
Af( 3 >(0) = - 



2 2d - 1 T 2 {d) 
■KaT(2d) ' 

2 2d - 1 r 2 (d) [l-#(d)] 



7ra 2 r(2d) 

2 2d r 2 (rf)y>(d) 

7ra 3 r(2d) ' 
6 2 2d r 2 (d) [<oip(d) - 6dip 2 (d) ~ dip^(d)] 
2 7ra 4 r(2d) 



4 2 2rf r 2 (d) [6^(d)+^(d)] 

[ ' ira 5 T(2d) 



where tfj^(z) = d n+1 logr(z)/d2 n+1 are the Polygamma functions. 
1.75 r 

1.5 \ 

1.25 ', 



(72) 




0.25 0.5 0.75 1 1.25 1.5 1.75 2 X/ ° 

x/<T 

FIG. 9. Comparison of best fit of Meixner distribution to truncated Levy distribution one by using the same 
a 2 and the kurtosis k (short dashed) and once using the same value and curvature at the origin (long-dashed). 
The parameters are a 1 = 0.280 and K = 12.7 as in the left-hand cumulative distribution in Fig. [jj The Meixner 
distribution with the same a 2 and k has parameters a = 2.666, d = 0.079, b — 0, the distribution with the same 
value and curvature at the origin has a = 0.6145, d = 1.059, 6 = 0. The very large cr-regime, however, is not fitted 
well as can be seen in the cumulative distributions which reach out to x of the order of 10ct in Fig. n. 



D. Other Non-Gaussian Distributions 

As mentioned in the beginning, some processes are also fitted quite well be generalized hyperbolic 
distributions They read 



and 



H G {x) 



( 7 2 -/3 2 ) V V' 
7*- 1 /2<5A V2^ 



2 i J2' 



A/2-1/4 K X -\/2 ('JV / S 2 +X 2 ) 



(73) 



G(p) 



6 X ( 7 2 - P 2 



^A/2 KxUj-f-ifl + ipf 



the latter defining another Hamiltonian 



(74) 
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H G ( P ) EE- log G( P ). 



(75) 



Introducing the variable £ = 8\/ 7 2 — (3 2 , this can be expanded in powers of p as in Eq. (pq), yielding the 
first two cumulants: 



ci = /3 



5 2 #i+a(C) 



C2 



Kx(C) 



Using the identity 



2v 

K v+1 {z) - K v -x(z) = —K v {z), 



the latter equation can be expressed entirely in terms of 

#i+a(C) 



r = r(C) 



Kx{0 



(76) 
(77) 

(78) 
(79) 



C2 



#i+a(C) 



*a(0 



b 2 fi 2 fi 4 

-'' + ^-[C + 2(l + A)r-Cr 2 ] 



c 



c 3 



(80) 



Usually, the asymmetry of the distribution is small implying that c\ and thus f3 is small. It is then useful 
to introduce the symmetric variance 



S 2 r/C, 



and write 

ci = f3a s , c 2 = er 2 = a 2 + (3 
The cumulants c 3 and C4 are most compactly written as 
c 3 =P 



^ + 6(l + \) d -a 2 -3ot 



+ (3 3 {2 (2 + A) ^ + [4 (1 + A) (2 + A) - 2( 2 ] ^a 2 - 6 (1 + A) ^ a 4 s + 2a s 6 | 



and 



, 3S 4 65 2 , N , n 4 
c 4 = «a 4 = — + — (1 + A) a 2 - 3d 4 



+ 6/3 2 |2(2 + A)^+[4(l + A) (2 + A) - 2( 2 ] ^ - 6 (1 + A) ^ o* a + 2a s 6 
+ /3 4 |[4(2 + A)(3 + A)-C 2 ]^+[4(l + A)(2 + A)(3 + A)-2(5 + 4A)C 2 ]^ 2 
- 2 [(1 + A) (11 + 7A) - 2C 2 ] ^cr 4 + 12 (1 + A) ^ - 3of| . 

The first term in C4 is equal to er^ times the kurtosis of the symmetric distribution 

3<5 4 6S 2 . n 
+ -r-s (1 + A) - 3, 



CM C 2 ^ 2 



(81) 



(82) 



(83) 



(84) 



(85) 
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Inserting here o 2 from (j3l|), we find 



«^^)+( 1 + A )^)- 3 - (86) 

Since all Bessel functions K v {z) have the same large- z behavior K v {z) — > y/n/2ze~ z and the small- z 
behavior K v (z) — > r(^)/2(z/2) !/ , the kurtosis starts out at 3/A for ( = and decreases monotonously to 
for £ — > oo. Thus a high kurtosis can be reached only with a small parameter A. 

The first term in C3 is (3K s aj, and the first two terms in C4 are K s aj + 6 (C3//3 — k s ct^) . For a symmetric 
distribution with certain variance a 2 and kurtosis n s we select some parameter A < 3/k s , and solve the 
Eq. ( p6| ) to find (. This is inserted into Eq. ( |8j] ) to determine 

,2 ^ 2 c (87) 



r(C) 

For larger kurtosis, the kurtosis is not an optimal parameter to determine generalized hyperbolic distri- 
butions. A better fit to the data is reached by reproducing correctly the size and shape of the near the 
peak and allow for some deviations in the tails of the distribution, on which the kurtosis depends quite 
sensitively. 

For distributions which are only slightly asymmetric, which is usually the case, it is sufficient to solve 
the above symmetric equations and determine the small parameter f3 approximately by the skew s = C3 /a 3 
from the first line in (^3|) as 



n s a s 



This approximation can be improved iteratively by reinserting (3 into the second equation in (|82]) to 
determine from the variance of the data a 2 an improved value of a 2 , then into the first two lines of 
Eq. (§3) to determine from the kurtosis k of the data an improved value of n s , and so on. 

E. Debye- Waller Factor for Non-Gaussian Fluctuations 

An important quantity for fluctuating Gaussian variables is the . It was introduced in solid state 
physics to describe the reduction of intensity of Bragg peaks due to the thermal fluctuations of the 
atomic positions. If u(x) is the atomic displacement field, the Gaussian approximation to the Debye- 
Waller factor e~ 2W is given by 

e- W = <e- V u(x) > = e- S *< l^'"*)] 2 )^ ( 89 ) 

This is a direct consequence of the relation for Gauss distributions 



(e Px ) 



dx —^ e - x2 / 2 ° 2 e Px = [dx e -V/2 e ^+P* = e ^P 2 /2 (90) 

J J 2ir K J 



which is a manifestation of Wick's rule 



whose right-hand side may be considered as the Debye- Waller factor for the Gaussian variable x. 
There exists a simple generalization of this relation to non- Gaussian distributions, which is 

P*\ = J dx~L { a T\x)e Px = j dx J ^. e -H( P ) e i P x+Px = p -H(iP) 



(91) 



e px ) = / dxL^' a '(x)e rx = dx ^r H We^ = e~ H ^ . (92) 
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F. Path Integral for Non-Gaussian Distribution 



Let us calculate the properties of the simplest process whose fluctuations are distributed according to 
any of the general non-Gaussian distributions. We consider the stochastic differential equation 



±{t) 



(93) 



where the noise variable rj(t) is distributed according an arbitrary distribution. The constant growth rate 
r x in (^3|) is uniquely defined only if the average of the noise variable vanishes: (r)(t)) = 0. The general 
distributions discussed above can have a nonzero average (x) = c\ which has to be subtracted from rj(t) 
to identify r x . The subsequent discussion will be simplest if we imagine r x to have replaced c\ in the 
above distributions, i.e., that the power series expansion of the Hamiltonian (fl9j) is replaces as follows: 



H{p) -> H r ^ (p) = H(p) - H'{0)p + ir x p = H(p) 



ir x p 



ir x p- 



1 



-C2P 



■ 1 3 

3! 3F 



1 



7C4P 



. 1 5 



(94) 



(95) 



2 * r 3!"°" 4! 

Thus we may simply work with the original expansion ( [49| ) and replace, at the end, 

ci -> r x . 

The stochastic differential equation ( |93|) can be assumed to read simply 

±(t) = v(t), (96) 

With the ultimate replacement ( |95| ) in mind, the probability distribution of the endpoints for the paths 
starting at a certain initial point is given by a path integral of the form 



P(x b t b \x a t a ) = I Vq I Vxexp 



dtH(r}(t)) 



S[x — 77], 



(97) 



with the initial condition x(t a ) — x a . The final point is, of course, x b = x{tb). 

The function H(r]) to be modified at the end by the replacement (j^) is the negative logarithm of 
the truncated Levy distribution or the generalized hyperbolic distribution or any other distribution. By 
analogy with the definition in momentum space (pi]) we shall define H (77) by 



e-^'si^ti), or e 
or for any distribution function D(x) generically by 

e -H(x) = 



-H{x) 



G(x), 



(98) 



(99) 



The correlation functions of the noise variable r)(t) in the path integral (|97|) are given by a straight- 
forward functional generalization of formulas (|53|). For this purpose, we express the noise distribution 



Pfo] = exp ~£dtH(r,(t)) 



P[rj\ 



Vr] 



]) as a Fourier path integral 
Vp 



2n 



exp 



[<p(t)»;(t) - H( P (t))] , 



(100) 



and note that the correlation functions can be obtained from the functional derivatives 



After n partial integrations, this becomes 



5p{h 



5 i f tb dtp(t)rj<t) 

P. J ta 



5p(t r . 



-J^dtH(p(t)) 



6 P (h 



8p(ti) Sp(t n ) 



S - ( tb dtH(p(t)) 

, f> <J t a 



5p(t 



p(t)=0 



6 - f tb dtH(p(t)) 

e Jt » 



(101) 
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f b dtH{p{t)) . 

J ta IT 



By expanding the exponential e 



in a power series we find immediately 



dtH(rj(t)) 



to. 



= Cl, 



Vrj 77(^1)77(^2) exp 
-h) + 4, 

v(ti)v(h)r](t3) exp 



dtH(r)(t)) 



dtH(r)(t)) 



= C 3 5(7j 1 -t 2 )5(ti-t 3 ) 

+ c 2 ci [*(*i-*a) + S(* 2 -t3)+<Kti-t3)] + c ?> 
(7/(^)^2)^3)^4)) = Z- 1 / P7y77(7; 1 )7 ? (i2)7 ? (7; 3 )7 ? (7i4)exp [- [" dtH(r)(t)) 



(102) 



(103) 



(104) 



= ci«(ii-t 2 )J(ti-ta)5(t 1 -i 4 ) 

+ C3C1 [S(ti—t2)S(ti — t 3 ) + 3cyclic perms] (105) 

+ 4 [*(ti-t2)*(t3-*4) + *(*l-*s)*(*a-*4) + «5(*l-*4)*(*2-*3)] 

+ c 2 Ci [<5(f 1 -7j 2 ) + 5 pair terms] + c\, 



where 



V-q exp 



dtH(r](t)) 



The higher correla tion fu nctio ns ar e ob vious generalizations of 
right-hand side of ( |l03| )-( |l05| ) Eq. ( |105| ) are distinguishable by t 



(106) 

The different contributions on the 
leir connectedness structure. 



An important property of the probability (|97|) is that it satisfies the semigroup property 

poo 

P{x c t c \x a t a ) = / dx b P(x c t c \x b t b )P(x b t b \x a t a ). 



(107) 



In Fig. 10 we show that this property is satisfied by experimental asset distributions reasonably well, 
except in the low-p robability tails. The discrepancy manifests itself also at another place: We shall see 
in Subsection III G that the solution of the path integral has a kurtosis decreasing inversely proportional 
to the time. The data in Fig. |l0|, however, show only an slower inverse square-root falloff. This can be 
accounted for in the theory by including fluctuations of the width a, which are certainly present as was 
illustrated before in Figs. |^ and [|. For calculations of this type with Gaussian distributions see Ref. p4| . 
If the semigroup property was satisfied perfectly, the Levy parameter a would be time independent as we 
can see from Eq. ( |4^ ) with a 2 oc (t b — t a ) and n oc l/(t b — t a ). With the slower falloff of n cx 1/ \/Jt b — t a ), 
however, a decreases like \j ^f(t b — t a ). 




k ps const, x 1/(4 — t a ) 



0.5 



see extra fig6b.pdf 
1G 



{tb - *a)/5 min 



FIG. 10. Left-hand side: Cumulative distributions obtained from repeated convolution integrals of the 15-min 
distribution (from Ref. |^|). Apart from the far ends of the tails, the semigroup property (107) is reasonably well 
satisfied. Right-hand side: Falloff of kurtosis is slower than l/(t& — t a ) expected from convolution property . 



G. Fokker-Planck-Type Equation 



The ^-functional may be represented by a Fourier integral leading to 



P(x b t b \x a t a ) 



dt 



ip{t)x{t) - ip{t)r)(t) - H( V (t)) 



(108) 



Integrating out the noise variable r)(t) amounts to performing the inverse Fourier transform (|2^) at each 
instant of time and we obtain 



P(x b t b \x a t a ) 



Vp 
2^ 



exp 



dt Mt)x(t) - H(p(t))} 



(109) 



Integrating over all x(t) with fixed end points enforces a constant momentum along the path, and we 
remain with a single integral 



(110) 



(111) 



/dp 
— exp [-(t b - t a )H(p) + ip(x b - x a )) ■ 

The Fourier integral can now be performed and we obtain, for a truncated Levy distribution: 

P(x b t b \x a t a ) = L^u^_ t \ix b — x a ). 

The result is therefore a truncated Levy distribution of increasing width. All expansion coefficients c„ of 
H{p) in Eq. ( |35| ) receive the same factor t b — t a , which has the consequence that the kurtosis n — c^/c^ 
decreases inversely proportional to t b —t a . The distribution becomes increasingly Gaussian with increasing 
time, as a manifestation of the central limiting theorem of statistical mechanics. This is in contrast to 
the pure Levy distribution which has no finite width and therefore maintains its power falloff at large 
distances. 

From the Fourier representation ( 110 ) it is easy to prove that this probability satisfies a Fokker-Planck- 
type equation 



d t P{x b t b \x a ta) = -H(-id x )P(x b t b \x a t a ). 



(112) 



Indeed, the general solution ip(x, t) of this differential equation with the initial condition ip(x, 0) is given 
by the path integral generalizing 



tp{x ) t) = / Vrj exp 



dtH(r)(t)) 



^[x- dt'r)(t') 



(113) 



To verify that this satisfies indeed the Fokker-Planck-type equation (112) we consider ip( x , t) a t a slightly 
later time t + e and expand 



i/j(x, t + e) 



T>rj exp 
T>r\ exp 
ij}'(x- 



dtH(rj(t)) 



dtH{r)(t)) 

ct+e 



If \^X — 
ijj ( X 



dt'ri(t') - 
dt'r)(t r 



t+e 



dt'rj(t') 



3!^ ' X 



dt'r){t') 
dt'rj(t' 
dt'rj(t' 



dt'r](t') 



t+e 



dtidt 2 r](ti)r)(t 2 ) 



(114) 



t+e 



dt 1 dt 2 dt 3 r n(t 1 )rj{t 2 )ri(t 3 ) 



+ ^ (4) (x - £ dt'ri(t')\ J dhdt 2 dt 3 dU v(hMhMh)v(ti) 
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Inserting here correlation functions (102)-(105) we obtain 



ip(x, t + e) 



T>r\ exp 



dtH(r)(t)) 



-e Cl d x + (ec 2 + e 2 ci) -d 2 x - (ec 3 + 3e 2 c 2 ci) -<9 3 



(115) 



+ (ec 4 + e 2 4c 3 ci + e 2 3cj + e 3 c 2 c 2 + e 4 c 2 )— <^ + 



■0 I a; 



dt'r)(t') 



In the limit e — > 0, only the linear terms in e contribute, which are all due to the connected parts of 
the correlation functions of rj(t). The differential operators in the brackets can now be pulled out of the 
integral and we find the differential equation 



d t t()(x,t) = 



- Cl d x + c 2 i<9 2 - c 3 ^d 3 x + c 4 ^d x 



ip (x, t) . 



(116) 



We now replace c\ — > r x and express using ( p^ ) the differential operators in brackets as Hamiltonian 
operator —H ra ,(~id x ). This leads to the Schrodinger-like equation 



d t ip(x,t) = -H r „(-id x )ij)(x,t) 



(117) 



Note that due to the many derivatives in H(id x ), this equation is in general non-local. 

By a similar procedure as in the derivation of Eq. (116) it is possible to derive a generalization of Ito's 
rule (|^) to functions of noise variable with non-Gaussian distributions. Thus we expand f(x(t + e)) as 
in (Q) where x(t) satisfies now the stochastic differential equation x(t) = r](t) with a nonzero expectation 
value (r](t)} = c\. In contrast to the previous evaluation of (^) with Gaussian noise, which had to be 
carried out only up to second order in the noise variable, we must now keep all orders. Evaluating the 
noise averages of the multiple integrals on the right-hand side using the correlation functions (102)— ( |l05| ), 
we find the time derivative of the expectation value of an arbitrary function of the fluctuating variable 
x(t) 



(f(x(t+e)))=(f(x(t))) 



-cid x 



{f( X {t)))ec 1 + -{I"{x{t))){ec 2+ e 2 c\) 
|</ (3) (*(*))>( 



ec 3 + e 2 c 2 ci + e 3 c?) 



(f(x(t)))+0(e 2 ). 



After the replacement c\ — > r x the function f(x(t)) obeys therefore the following equation: 

(f(x(t)))=-H rx (id x )(f(x(t))). 
Taking out the lowest-derivative term this takes a form 

(f(x(t))) = (d x f(x(t)))(x(t))-H rx (id x )(f(x(t))). 



(118) 



(119) 



(120) 



In postpoint time slicing, this may be viewed as the expectation value of the stochastic differential 
equation 



f(x(t)) = d x f{x{t))x{t) - H rx (id x )f(x(t)). 

This is the direct generalization of Ito's rule (||). 

For an exponential function f(x) = e Px , this becomes 



dt 



> P "M = [x(t)-H rx (iP)]e Px W. 



As a consequence of this equation for P = 1, the rate rs with which a stock price S(t) 
on the average according to formula (fD) is now related to r x by 



c(t) 



(121) 



(122) 



grows 
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r s =r x - H(i) = r x - [H{i) - iH'(0)] = -H rm (i), (123) 

which replaces the simple Ito relation r$ — r x + a 2 /2 in Eq. (|Io|). Recall the definition H(p) = H(p) — 
H'(0)p in Eq. (94). The corresponding generalization of the left-hand part of Eq. (||) reads 



| = x(t) - H(i) = x(t) - [H{i) - iH'(0)} = x{t) -r x - H rx (*). (124) 



The forward price of a stock must therefore be calculated with the generalization of formula (|llj) : 

(S(t)) = S(0)e rst = S(0)(e r * t+ f« dt ' vit ' } ) = S(0) e - ff - Wt = S (0)e^ H ^- lH ' . (125) 



Note that e may derive the differential equation of an arbitrary function f(x(t)) in Eq. ( |119| ) from a 
simple mnemonic rule, expanding sloppily 

f(x(t + dt)) = f(x(t) + xdt) = f(x(t)) + f'(x(t))x(t)dt + ±f"(x(t))x 2 (t)dt 2 

+ ^f {3) (x(t))x 3 (t)dt 3 + ... , (126) 

and replacing 

{x(t))dt ->■ adt, (x 2 (t))dt 2 -» c 2 dt, (x 3 (t))dt 2 -> c 3 dt, .... (127) 



IV. CONCLUSION 

The new stochastic calculus developed in this paper should be useful for estimating financial risks of a 
variety of investments. In particular, it will help in developing simple methods to estimating fair option 
prices for stocks with non-Gaussian fluctuations. More details can be found in the textbook Ref. pgfl. 
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